/*
 * OdeSolver.h
 *
 *  Created on: 28 Jul 2011
 *      Author: Allan
 */

#ifndef ODESOLVER_H_
#define ODESOLVER_H_

// C++ includes
#include <memory>
using namespace std;

// Eigen includes
#include <Eigen/Dense>
using namespace Eigen;

class OdeProblem
{
public:
	virtual void Function(double t, const VectorXd& y, VectorXd& f) = 0;

	virtual void Jacobian(double t, const VectorXd& y, MatrixXd& dfdy, VectorXd& dfdt) {};
};

enum OdeAlgorithm
{
	RungeKutta2,
	RungeKutta4,
	RungeKuttaFehlberg45,
	RungeKuttaCashKarp45,
	RungeKuttaPrinceDormand8,
	RungeKuttaImplicitGaussian1,
	RungeKuttaImplicitGaussian2,
	RungeKuttaImplicitGaussian4,
	BulirschStoerImplicit,
	MultistepAdams,
	MultistepBDF
};

class OdeSolverImp;

class OdeSolver
{
public:
	OdeSolver();
	
	OdeSolver(unsigned dimension, OdeAlgorithm algorithm = RungeKutta4);
	
	~OdeSolver();
	
	void Initialiaze(unsigned dimension, OdeAlgorithm algorithm = RungeKutta4);
	
	void SetControlStandard(double epsAbs, double epsRel, double yFactor, double dydtFactor);
	
	void SetControlSolution(double epsAbs, double epsRel);
	
	void SetControlSolutionDerivative(double epsAbs, double epsRel);
	
	bool Iterate(OdeProblem& problem, double& t, double tEnd, double& dt, VectorXd& y);
	
	bool Solve(OdeProblem& problem, double& t, double tEnd, double& dt, VectorXd& y);

public:
	std::unique_ptr<OdeSolverImp> imp;
};

#endif /* ODESOLVER_H_ */

